
VORTICITY INTERACTION EFFECTS 
ON BLUNT BODIES 


E, Clay Anderson and David C. Wilcox 

Prepared by 
DCW INDUSTRIES. 

Sherman Oaks, Calif. 91423 
for Langley Research Center 



V: 



NATIONAL AERONAUTICS AND SPACE ADMINISTRATION • WASHINGTON, D. C. • JANUARY 1977 


TECH LIBRARY KAFB, NM 



I 


1. Report No. 

NASA CR-2778 

4. Title and Subtitle 


2. Government Accession Mo. 


VORTICITY INTERACTION EFFECTS ON 
BLUNT BODIES 


7. Author(s) 

E. Clay Anderson 
David C. Wilcox 


9. Performing Organization Name and Address 

DCW Industries 

13535 Ventura Boulevard, Suite 207 
Sherman Oaks, CA 91423 


12. Sponsoring Agency Name and Address 

NASA Langley Research Center 
Hampton, VA 23365 

15. Supplementary Notes 

Langley technical monitor: Kenneth Sutton 

Final report. 

16. Abstract 


TECH LIBRARY KAFB, NM 


Y 8. } Performing Organization Report No. 

Dcw-R-08-01 


3. Recipient's Catalog No. 


□ □1*1383 


5. Report Date 

January 1977 


6. Performing Organization Code 

_ DW-08 


10. Work Unit No. 


11. Contract or Grant No. 


nasi-1397^ 

13. Type of Report and Period Covered 

Final Report - 
-6/1/75-6/1/7^ 

14. Sponsoring Agency Code 


A two-part study of hypersonic viscous shock layers has been conducted. 
In the first part of the study, numerical solutions of the viscous-shock- 
layer equations governing laminar and turbulent flows of a perfect gas 
and radiating and non-radiating mixtures of perfect gases in chemical 
equilibrium are presented for hypersonic flow over spherically blunted 
cones and hyperboloids. Turbulent properties are described in terms of 
the classical mixing length. Results are compared with boundary-layer 
and inviscid flowfield solutions; agreement with inviscid flowfield data 
is satisfactory. Agreement with boundary-layer solutions is good except 
in regions of strong vorticity interaction; in these flow regions, the 
viscous-shock-layer solutions appear .to be more satisfactory than the 
boundary-layer solutions. The second part focuses upon preparation for 
future applications; boundary conditions suitable for hypersonic viscous 
shock layers are devised for an advanced turbulence theory. 


17. Key Words (Suggested by Author(s) ) 


18. Distribution Statement 

Viscous Shock Layer; 

Radiation 

Unclassified - Unlimited 

Boundary Layer; 

Gasdynamics 


Turbulence Model; 
Streamline Swallowing; 
Vorticity Interaction; 


Subject Category 02 


19. Security Gassif. (of this report) 

20. Security Classif. (of this page) 

21. No. of Pages 

22. Price* 

Unclassified 

Unclassified 

70 

$4 . 25 


For sale by the National Technical Information Service, Springfield, Virginia 22151 



CONTENTS 


SECTION PAGE 

NOTATION v 

1 INTRODUCTION 1 

2 MIXING LENGTH APPLICATIONS 5 

2.1 Governing Equations 5 

2.1.1 Boundary Conditions 8 

2.1.2 Radiative Transport 9 

2.1.3 Thermodynamic and Transport 

Properties 9 

2.2 Eddy-Viscosity Approximations 10 

2.2.1 Inner-Eddy-Viscosity 

Approximation 10 

2.2.2 Outer-Eddy-Viscosity 

Approximation 11 

2.3 Solution Method 13 

2.4 Results and Discussion 17 

2.4.1 Perfect Gas Solutions with 

Strong Vorticity Interaction 17 

2.4.2 Solutions for Radiating and 
Nonradiating Gas Mixtures in 

Chemical Equilibrium 24 

3 DEVELOPMENT OP AN ADVANCED TURBULENCE MODEL ... 35 

3.1 The Model Equations 36 

3.2 Surface Boundary Conditions 37 

3.3 Turbulent-Nonturbulent Interface 

Structure 43 

3-3-1 Interface Structure for a 

Conventional Boundary Layer 43 

3-3.2 Interface Structure for 

Stagnation Point Plow 45 

3-3.3 Tentative Modification to the 

Model Equations 48 



CONTENTS (concluded) 


SECTION PAGE 

3 DEVELOPMENT OF AN ADVANCED TURBULENCE 

MODEL (continued) 

3.4 Computing-Time Reduction 53 

H CONCLUDING REMARKS 58 


REFERENCES 


59 



NOTATION 


SYMBOL 

A + 


C I 


P 


P 


c’ 

e 

H 

H, 


h 

h 

h. 


A 




K 

k 

k + 

£ 

£ 

M* 

M* 

N 


DEFINITION 

Damping factor [Equations (15) and (l6)J 
Mass fraction of species i, p^/p 
Mass fraction of element £ 

N 

Frozen specific heat of mixture, X C. C . 

i-1 p ’ 

Specific heat of species i, C* ^/C* 

Inviscid flow velocity gradient 

Turbulent mixing energy 

u 2 

Defined quantity, h + — 

v 2 

Total enthalpy, H + 
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1=1 1 1 
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VORTICITY INTERACTION EFFECTS ON BLUNT BODIES 


E. Clay Anderson* and David C. Wilcox"^" 

DCW Industries 

1 . INTRODUCTION 

Vorticity interaction between the boundary layer and the outer 
inviscid entropy layer is not significant for high Reynolds 
number laminar flows over blunt bodies having low angle after- 
bodies except at distances far downstream. However, for large 
angle afterbodies such as those being considered for planetary 
entry probes, 1 vorticity interaction is strong in the region 
near the nose of the body if the flow is turbulent. For free- 
stream conditions which require an analysis including mass 
injection at the surface, vorticity interaction can be a sig- 
nificant influence for laminar flows if the mass injection rate 
is sufficiently large. As a result of the high speeds asso- 
ciated with typical planetary entry trajectories, radiation heat 
transfer effects must also be included in the flowfield analysis 

The use of higher order boundary-layer theory for the analysis 
of flowfields with a high degree of coupling between the bound- 
ary layer (BL) and the inviscid flow requires a complex intera- 

tive solution procedure. Because of the difficulties experi- 

2 R 

enced in the application of higher order BL theory, Davis 5 

has developed a numerical technique for solving the viscous- 

shock-layer (VSL) equations governing laminar flow of a perfect 

gas and for binary mixtures with finite rate chemical reactions. 
2 | 

Moss has developed VSL solutions for multicomponent gas mix- 

R 

tures for both equilibrium and finite rate chemistry. Moss has 

* ~ 

Consultant 

+ 

Owner 


extended the equilibrium chemistry analysis to include the 
effects of radiation heat transfer. These YSL analyses have 
been found accurate for flows with strong vorticity interaction 
and/or radiation heat transfer. 

g 

In recent publications by Eaton and Larson 0 and Anderson and 
7 

Moss, numerical solutions of the YSL equations governing tur- 
bulent flow of a perfect gas have been presented. Eaton and 
Larson considered the thin VSL equations and presented 
solutions for turbulent flow over slender cones. Their 
results were compared with experimental data and showed good 
agreement. Anderson and Moss considered turbulent flow over 
blunt axially symmetric bodies and used the full VSL equa- 
tions. These solutions were compared with first-order tur- 
bulent BL theory. The results obtained with this turbulent 
BL analysis were essentially identical to first-order turbu- 
lent BL results in the region where vorticity interaction is 
not significant. For the downstream region, the expected 
differences were obtained. 

g 

In a later investigation, Anderson and Moss extended the 
turbulent VSL analysis to reacting gas mixtures in chemical 
equilibrium and considered both radiating and nonradiating 
shock layers. These turbulent VSL solutions were compared 
with first-order BL solutions and first-order BL theory with 
corrections for streamline swallowing (BLSSW). A nonradiating 
shock-layer solution corresponding to a typical trajectory 
point for a Venusian entry was considered. A coldwall bound- 
ary condition was assumed, and the results of the VSL analy- 
sis were compared with integral BL and BLSSW analyses. The 
two BL solutions and the VSL solution showed excellent agree- 
ment in the laminar flow region. For the turbulent flow 
region, both the BLSSW and the VSL solution showed a large 
increase in the heat-transfer rate when compared with the BL 
solution. For the turbulent flow region, the VSL heat-transfer 
2 


rates were 15$ to 30 $ lower than the solution obtained with 
the Integral BLSSW analysis, but the general trend In the 
heat-transfer rate distribution was the same for both methods 
of analysis. 

For the more complex case of a radiating shock layer with mass 
injection at the surface, a typical trajectory point For a 
Venusian entry was analyzed. These VSL results were compared 
with an integral matrix BLSSW analysis. For this case, the 
laminar flow region solutions differed by 10$ to 15$ in con- 
vective heating rate predictions, and for the turbulent flow 
region, the two methods of analysis showed opposite trends in 
the convective-heating-rate distribution. The BLSSW analysis 
predicted a maximum convective heating rate at the transition 
point. The VSL analysis predicted the maximum convective 
heating rate at a downstream location. For the predominantly 
inviscid outer flow region of the shock layer, the VSL solu- 
tion was found to be in satisfactory agreement with the invis- 
cid flowfield solution. The prediction of radiation heating 
at the surface was in satisfactory agreement with the BLSSW 
- analysis . 

The opposite trends in the calculated convective heating-rate 
distributions appeared to be the result of an increase in vor- 
ticity interaction caused by mass injection and the higher 
temperature at the surface boundary. This behavior was not 
obtained for the nonradiating case with no injection and cold- 
wall boundary conditions. The coldwall boundary condition 
results in a thinner BL thickness and diminishes the influence 
of vorticity interaction. Both higher surface temperature and 
mass injection increase vorticity interaction. 

The present investigation is directed primarily toward deter- 
mining the behavior of BLSSW and VSL solutions in regions of 
strong vorticity interaction. The study consists of two parts. 



In the first part of the study, a perfect gas is considered, 
and the eddy viscosity modeling techniques are identical in 
both the VSL and BLSSW analyses. Representative VSL profile 
data are compared with two different inviscid flowfield solu- 
tions and BL and BLSSW solutions. YSL solutions are also pre- 
sented for both radiating and nonradiating mixtures of perfect 
gases in chemical equilibrium, but comparative BL and BLSSW 
data are not available for the radiating shock-layer calcu- 
lations. The second part of the study focuses upon develop- 
ment of an advanced turbulence model in preparation for 
future viscous shock layer applications. 


4 



2. MIXING LENGTH APPLICATIONS 


In this section, hypersonic viscous shock layer solutions are 
presented. The equations of motion appropriate for shock 
layers are stated, including boundary conditions and assumed 
thermodynamic properties. The mixing-length closure approxi- 
mations are then discussed. The solution method is described 
next. Results of numerical computations are then presented 
for both perfect-gas and equilibrium radiating gas mixture 
applications . 

2.1 GOVERNING EQUATIONS 

The equations of motion for reacting gas mixtures in chemical 

q 

equilibrium are presented by Bird, Stewart, and Lightfoot. 

The formulation of these equations in body-oriented coordinates 

appropriate for viscous-shock-layer analysis of laminar flow 

4 5 

of radiating and nonradiating gases is presented by Moss. * 

For turbulent flow, the viscous-shock-layer equations are 
derived using methods analogous to those presented by Dorrance 10 
for the turbulent-boundary-layer equations and are expressed 
in nondimensional form for the coordinate system shown in 
Figure 1 as : 


Continuity : 


•sfg-C (r + n cos 0 ) puj + t^-[ (1 + nic ) (r + n cos 0 ) pvj = 0 (1) 


s- momentum: 


u 3u , 3u , uvk 
3n 1+ntc 


p \l+nic 3s 


-!_!£ = s 


1+nK 3s 


_ ~2; 3 T /--i . +-\ PUK 1 

~ a |toL ,,(1+e 1 to-rFs?J 


, (2k, 3 cos 0 NT , + 

\l+nK r + n cos 0/[_^ e 


) _ _M^|> ( 2 ) 

; 3n l+n< u K J 
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n-momentum : 


. / u 3v . ,.3v u 2 k \ . 3p _ n 
p \l+n«: 8s 3n - 1+nic ) 3n 


(3) 
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(4a) 


(4b) 


(4c) 


(4d) 


and 


H E h + y 


(4e) 


Elemental continuity: 


pU ^ ♦ V- 3C ^ 


,l+n< 3s 


3n 


(l+n<) (r + n cos 9 )<J 


3n 


(1+nK) (r + n cos 9 ) (n t + e + Pr — N 


3C , 


where 


N pr \ x, Le ’ & N prjT ”Le,T/ 3n 


N M 


Z 6. -|c. 

1=1 1 


(5a) 


(5b) 
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State : 


p = pTR*/M*C* (6) 

P,oo 


2.1.1 Boundary Conditions 

The boundary conditions at the shock are calculated -by using the 
Ranklne- Hugoniot relations. At the wall, the no-slip and no- 
temperature-jump boundary conditions are used; consequently, 
u w = 0 . The wall temperature and mass Injection rate are either 
specified or calculated. For the calculated mass injection con- 
ditions, the ablation process is assumed to be quasi-steady and 
the wall temperature is the sublimation temperature of the abla- 
tor surface. With these assumptions, the expression for the 
coupled mass injection rate is 


m 



* 

c ,w 


(C i h i ) 




(7) 


For ablation injection, the elemental concentrations at the wall 
are governed by convection and diffusion as given by the equa- 
tion 



Precursor effects are neglected while the energy reradiated from 
the surface is included in the radiation transport calculations. 
The net radiative flux, q , can be represented as the differ- 
ence of two components 


q 


r 



q. 


(-) 


(9) 
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At the surface 


q ( ‘ + ) * = eb*T * 4 (10) 

H r,w w 

where e is the surface emissivity. 

The heat transferred to the wall due to conduction and diffusion 
is 

( 11 ) 


-q 


C jW 


52 Kh + 


N 


N 
£ 

Sc i=l 


h . 

i 


SC. 

i 

3n 


2.1.2 Radiative Transport 

The radiative flux, q , is calculated with the radiative trans- 
it 12 r 

port code RAD. 5 The RAD computer code has been incorporated 
in the present viscous-shock-layer computer code (HYVIS) and 
streamlined for computational efficiency. 

The RAD code accounts for the effects of nongray self-absorption. 
Molecular band, continuum, and atomic line transitions are 
included. A detailed frequency dependence of the absorption 
coefficients is used for integrating over the radiation fre- 
quency spectrum and the tangent slab approximation is used for 
integrating over physical space. 

2.1.3 Thermodynamic and Transport Properties 

The equilibrium composition is determined by a free energy mini- 
mization calculation as developed "in Reference 13- Thermodynamic 
properties for specific heat, enthalpy, and free energy and 
transport properties for viscosity and thermal conductivity are 

required for each species considered. Values for the thermo- 

14 15 16 

dynamic s and transport properties are obtained by using 

polynomial curve fits. The mixture viscosity is obtained by 

17 

using the semiempirical formula of Wilke. 
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2.2 


EDDY-VISCOSITY APPROXIMATIONS 


A two-layer eddy-viscosity model consisting of an inner law 

based upon Prandtl's mixing-length concept and the Clauser- 

Klebanoff expression (based on References 18 and 19) for the 

outer law is used in the present investigation. This model, 

20 

introduced by Cebeci, assumes that the inner law is appli- 
cable for the flow from the wall out to the location where the 
eddy viscosity given by the inner law is equal to that of the 
outer law. The outer law is then assumed applicable for the 
remainder of the viscous layer. It is noted that the eddy 
viscosity degenerates to approximately zero in the inviscid 
portion of the shock layer. The degeneracy is expressed in 

19 

terms of the normal intermittency factor given by Klebanoff. 
The expressions used in the present investigation are given 
in the following sections. 


2.2.1 Inner-Eddy-Viscosity Approximation 

Prandtl’s mixing-length concept is stated in nondimensional 
variables as 


et = 

i aV 8nl 


( 12 ) 


~ 21 
The mixing length, £, is evaluated by using van Driest' s pro- 
posal stated as 


where 



(13) 


( 1 ^) 


Here, k 1 is the von Karman constant, which is assumed to have a 
value of 0.4, and A + is a damping factor whose value is 26. 


10 



20 

Cebeci suggests that for flows with a pressure gradient, the 
damping factor be expressed as 


A + = 26(1- 11 . 8P + ) 


-is 


(15) 


and for flows with both a pressure gradient and mass injection. 


\+ - 


26 <(- £i[exp(11.8v+) - 1] + exp(11.8v+) 


-h 


% 


( 16 ) 


w 


where 


■ - 52 (lf) 


e P 2 u* 


(17) 


v 

,+ = _W 
w u 


( 18 ) 


and 
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w 
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(t) 

\ / w 


(19) 


2.2.2 Outer-Eddy-Viscosity Approximation 

For the outer region of the viscous layer the eddy viscosity is 
approximated by the Clauser-Klebanof f expression 


£ 


+ 

O 


k„ pu 6, Y . 
2 K e k ' l 


a 2 ii 


3 n 


where k ^ 


0.0168 and 6, and y. 

k 'i,n 


are defined as 


6 


k 



e ' 


( 20 ) 


( 21 ) 
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and 



( 22 ) 


20 

Equation (22) is Cebeci's approximation of the error-function 

19 

definition presented by Klebanoff. 


For equilibrium flow without radiation, the boundary-layer thick- 
ness <5 is assumed to be the value of n at the point where 

H f 

j~- = 0.995 (23) 

t ,co 

and is defined by linear interpolation in an array of local total 
enthalpies. This definition is approximately equivalent to the 
usual boundary-layer definition 

~ = 0.995 (24) 

e 

where u is the local value for the undisturbed inviscid flow 
e 

outside the boundary layer. 

The values of the parameters k-^ and k 2 in Equations (13) and 

(20) depend on the flow conditions being considered, as does 

the constant represented by the value 26 in Equations (15) and 

(16). The values given are used for convenience in developing 

the numerical method. A discussion of these parameters is 

22 

presented by Harris. 


For radiating gases, the loss of energy from the shock layer 
makes the total enthalpy definition unsatisfactory. For these 
cases, the boundary-layer thickness is assumed to be defined 
by an index of diffusion, conduction, and dissipation expressed 
as 


12 



f o 1 

latwh) +:| 

[ K 

, .1 cos 0 N 

) ( 4 ) 1 +<} ) 2 + ^3 ) 

| dn 

\ 1+nc 

r + n cos Q; 


+ 1 

( * 

, ,1 cos 0 N 

) ( <J> 1 +6 2 +( f , 3) 

dn 

1+nic 

r + n cos 0, 


0.95 (25) 


Boundary-layer thicknesses given by this expression are essen- 
tially the same as those given by the total enthalpy definition 
when applied to perfect gas flows, but both definitions over- 
predict the BL thickness in the vicinity of the stagnation 
streamline when an equilibrium gas is considered. 


For the present study, the turbulent Prandtl and Lewis numbers 
are assumed to be 0.9 and 1.0, respectively. 

2.3 SOLUTION METHOD 

2 

Davis presented a method for solving the viscous-shock-layer 

4 5 

equations for stagnation and downstream flow. Moss 5 applied 
this method of solution to reacting multicomponent mixtures. 
The present method of solution is identical to that of Refer- 
ences 2 and 4. Therefore, only an overview of the solution 
procedure is presented here. 


The numerical computation is simplified by normalizing most of 
the variables with their local shock values. The transformed 
independent and dependent variables are 


n = n/n g 

P = P/P 

5 = s 

P = p/p 

u = u/u s 

T = T/T 

v = v/v s 

H = H/H 

u = u/y s 

K = K/K 


( 26 ) 


13 



Since the normal coordinate, n, is normalized with respect to the 
local shock standoff distance, a constant number of finite-differ- 
ence grid points between the body and shock are used. The trans- 
formations relating the differential quantities are 


_ 8 _ = _ 8 _ ^s _ 8 _ 

8s 8£ “ n n 8n 
s 


( 27 a) 


where 


and 



_ 8 _ _ _ 1 _ _ 8 _ 

8n n 8r| 
s 


3 2 = JL_ 8 2 

8n 2 n 2 8n 2 
s 


( 27 b) 


(27c) 


(27d) 


The second-order partial differential equations are linearized 
and written in the standard form for a parabolic equation as 

+ a l3W + a 2 w + a 3 + “4 35 ■ 0 (28) 

where W represents tangential velocity for the s-momentum equa- 
tion, enthalpy for the energy equation, and elemental concen- 
trations for the elemental continuity equations. For the energy 
equation, the divergence of the radiative flux is included in 
the term. When the derivatives in Equation (28) are con- 
verted to finite-difference form by using Taylor’s series expan- 
sions, the resulting equations are of the following form: 


AW , + B W 

n m,n-l n n 


(29) 


14 



The subscript n denotes the grid points along a line normal 
to the body surface, and the subscript m denotes the grid 
stations along the body surface. Equation (29), along with 
the boundary conditions, constitute a system of the tridiagonal 
form and can be solved efficiently. 

A variable grid spacing is used in both the tangential and 
normal directions to the surface so that the grid spacing can 
be made small in the region of large gradients. The order of 
the truncation terms neglected are AC (first order accurate) 
and either Ap Ap ^ or ( An n - Ap ^ ) . 

The equations are solved at any body station m in the order 
shown in Figure 2. The governing equations are uncoupled and 
the values of the dependent variables are computed one at a 
time. Each of the second-order differential equations is 
individually integrated numerically by using the tridiagonal 
formalism [Equation (29)]. The global continuity equation is 
used to obtain both shock standoff distance and the v compo- 
nents of velocity. By integrating the global continuity equa- 
tion between the limits of 0 to 1, a quadratic equation for n_ 

o 

is obtained. For the v component of velocity at p , the global 
continuity equation is integrated with respect to p between the 
limits of 0 to p. The pressure, p, is determined at station m 
by integrating the normal momentum equation with respect to p 
between the limits of 1 to p. The equation of state is used 
to determine the density. The solution is iterated until con- 
vergence is achieved. The solution advances to the next body 
station, m+1, and uses the previous converged solution profiles 
as initial values for starting the solution at station m+1. 

This procedure is repeated until a solution pass is obtained. 

The first solution pass provides a first approximation to the 
flowfield solution because the following assumptions are used. 

15 




Figure 2. Flow design of local solution procedure. 









The thin shock-layer form of the n-momentum equation, X+nl< = 8n" 2 * * 5 
is used; the stagnation streamline solution is independent of 
downstream influence; the term dn /d£ is set to zero at each 
body station, and the shock angle a is assumed to be the same 
as the body angle 0. These approximations are then removed 
by global iteration- Two solution passes are generally suffi- 
cient. This solution procedure is programmed for the CDC 6600 
computer . 

2.4 RESULTS AND DISCUSSION 

Numerical solutions of the VSL equations governing laminar, 
transitional, and turbulent flows of a perfect gas are compared 
with inviscid, BL, and BLSSW solutions. Equivalent but less 
extensive comparisons are presented for nonradiating mixtures 
of perfect gases in chemical equilibrium. VSL solution data 
are also presented for a radiating mixture of perfect gases in 
chemical equilibrium but comparative data are not available 
for this case. 

2.4.1 Perfect Gas Solutions with Strong 
Vorticity Interaction 

VSL solutions for laminar and turbulent flow of a Mach 10 free- 
stream over a 40° half-angle spherically blunted cone are pre- 
sented. Representative profile data obtained using the VSL 
analysis are compared with inviscid, BL, and BLSSW profile 
data. The inviscid flowfield solutions were obtained using 

23 

the blunt-body method of characteristics (MOC) of Inouye, et al. 

24 

and the transient finite-difference (TFD) procedure of Sutton. 

The BL data were obtained using the analysis of Anderson and 
25 

Lewis. BLSSW data were obtained using the analysis of Mayne 

2 6 

and Dyer (this solution was provided by Dr. Arloe W. Mayne, Jr. 

of ARO, Inc.), and the analysis of Reference 25 was modified to 

account for variable entropy edge conditions during the present 
investigation. 
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Surface pressure distributions corresponding to the VSL and 
MOC solutions are shown in Figure 3. Excellent agreement is 
obtained except in the vicinity of the sphere-cone tangency 
point where differences are ±10$. These differences have little 
influence upon the downstream flow and are not excessive in the 
region near the tangency point. The shock shapes corresponding 
to the VSL and MOC solutions are essentially identical and are 
not presented. 

Figure 4 shows heating-rate distributions corresponding to 
(1) the BL analysis of Reference 25, (2) the BLSSW analysis of 
Reference 26, (3) VSL solutions corresponding to the assumptions 
that the BL edge is located as defined by either Equation (23) 
or (25) and at the bow shock, and (4) laminar VSL and BLSSW 
(Reference 25) solutions. For turbulent flow, instantaneous 
transition at s=0.8 has been assumed, and the damping factor, 
A + , in the van Driest eddy viscosity law is assumed to be 26 in 
all calculations. 

The turbulent VSL solution with the BL edge location defined by 

either Equation (23) or (25) and using Equation (22) for y. 

a s n 

is in excellent agreement with the BL analysis of Reference 25 
for s£2.5* In the region of strong vorticity interaction, the 
heating rates predicted by the VSL solution are as much as 40% 
less than the heating rates predicted by the two BLSSW solu- 
tions. In the downstream region where the inviscid entropy 
layer has been swallowed by the BL, the VSL and BLSSW solutions 
approach an equivalent cone solution, as they should. 

The turbulent VSL solution is strongly influenced by the defi- 
nition of the BL edge location. To demonstrate the maximum 
influence of the BL edge definition on the VSL solution, the 
BL edge was assumed to be located at the bow shock and a unit 

value for y. was assumed. These assumptions result in maximum 
i , n 
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Figure 3- Surface pressure distributions. 
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eddy viscosities using the present mixing length turbulence 
model. This solution shows the expected behavior, but the peak 
heating rates predicted by the VSL solution remains substan- 
tially less than that given by the BLSSW solutions. The lami- 

25 

nar VSL and BLSSW solutions which are not directly influenced 
by the BL and the displacement thicknesses show excellent agree- 
ment (±3$) for the heating-rate distribution. The data presented 
are the VSL solution. 

The apparent contradiction for the agreement between laminar 
and turbulent VSL and BLSSW solutions is the result of the BL 
and the displacement thickness distributions obtained using the 
BLSSW analyses. Both BLSSW solutions showed a linear increase 
in the BL and the displacement thickness distributions up to 
the point where the inviscid entropy layer is swallowed by the 
BL. Downstream of this point, these thickness distributions 
showed a gradual decrease in magnitude until the equivalent cone 
solution was approached. Since these definitions determine the 
scale for turbulence, the turbulent-heat-rate distributions 
predicted by the BLSSW solutions show the same behavior. 

The rapid growth of the BL and displacement thickness in the 

region where longitudinal entropy gradients were significant 

was obtained for both laminar and turbulent flow when the vari- 

27 

able entropy corrections presented by Blottner were used in 
the analysis of Reference 25. An approximate streamline swallow- 
ing analysis which corrected only the edge conditions was also 

2 8 

used. This approach has been used by Price and Harris and was 
found to give a more satisfactory distribution for the BL and 
the displacement thicknesses. The heat-transfer-rate distribu- 
tion obtained with this approximate analysis was essentially 

27 

the same as when the full modifications presented by Blottner 

were used. A more accurate BL solution would require higher- 

2Q 

order theory (see van Dyke *). 
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Figure 5 shows velocity profile data at four locations. The VSL 
solutions are compared with the inviscid flowfield solutions of 
References 23 and 24, the BLSSW solution of Reference 26, and the 
BL solution of Reference 25. At s=2.1, the VSL profile data in 
the inner viscous region show good agreement with both BL and 
BLSSW solutions for turbulent flow. However, as is evident from 
the different heating rates predicted at this station, the temper- 
ature profiles corresponding to the BLSSW solution showed pro- 
nounced differences from that of the VSL solution. The tempera- 
ture profiles corresponding to the BL solution show good agreement 
with the VSL solution as indicated by the heating-rate prediction. 
These data are not presented. 

In the predominantly inviscid outer flow region, the VSL pro- 
file data show satisfactory agreement (±5%) with both inviscid 
flowfield solutions. The VSL and TFD solutions do not show the 
distinct sharp peak in the velocity profile that is a charac- 
teristic of spherically blunted cone flowfields as shown by 
the MOC results. Since TFD inviscid flowfield solutions gen- 
erally use less than 15 interior grid points across the shock 
layer, the differences between the TFD and MOC solutions are 
probably the result of inadequate resolution. The VSL solu- 
tions were obtained using 150 interior grid points with 
approximately 80 points within the predominantly viscous inner 
region and the remaining 70 points in the predominantly invis- 
cid outer region. Since velocity and temperature gradients 
are relatively large in the inviscid entropy layer, the differ- 
ences between the MOC and VSL solutions appear to be the result 
of viscous and displacement effects. 

Profile data at s = 3-2 and 6.7 show that the agreement between 
the VSL and BLSSW solution becomes progressively worse as vor- 
ticity increases at the BL edge. After the inviscid entropy 
layer is swallowed by the viscous layer, the agreement between 
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Figure 5. Comparisons of velocity profiles. 
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the two methods of solution is satisfactory as shown at s = 21. 
The profile data at s = 6.7 show that the inviscid entropy 
layer has been swallowed by the BL using the BLSSW analysis. 

The VSL profile data show a distinct "inviscid entropy layer," 
and these data are in satisfactory agreement with both invis- 
cid flowfield solutions. 

The solutions presented for the 40° half-angle spherically 
blunted cone indicate that the VSL solution procedure is satis- 
factory for high Reynolds number laminar and turbulent flows. 

The VSL solutions using algebraic mixing length turbulence 
models in conjunction with the definition expressed by either 
Equation (23) or Equation (25) for the BL edge location show 
acceptable agreement with BL solutions in the nose region of 
blunt bodies and with BLSSW solutions in the far downstream 
flow. The solutions appear to be more satisfactory than BLSSW 
solutions for flows with strong vorticity interaction. However, 
the accuracy of the present solution procedure and the appli- 
cability of the algebraic mixing length turbulence modeling 
technique in the presence of strong vorticity interaction must 
be determined by experimental verification. Turbulence models 
which are independent of conventional boundary-layer thickness 
definitions are needed for the present method of analysis. 

2.4.2 Solutions for Radiating and Nonradiating 
Gas Mixtures in Chemical Equilibrium 

A VSL solution for the flowfield over a 60° half-angle spheric- 
ally blunted cone for a typical Venusian entry trajectory point 
is presented in Figures 6 and 7- The VSL solution is compared 

with BL and BLSSW solutions obtained using the integral analy- 
30 

sis of Edquist (this solution was provided by C.T. Edquist, 
Martin-Marietta Corp. , Denver Division). 

The inviscid flowfield solution used to specify edge conditions 
for the BL and BLSSW solutions was determined using a single 


24 



60° SPHERE CONE 
U# = lO.Okm/s 

3 

= 0.01 kg/m 



I i i i i 

0 0.5 1.0 1.5 2.0 

BODY RADIUS, r 


Figure 6. Surface pressure distributions, 
ro 

VJl 


NONRADIATING EQUILIBRIUM GAS (0.95 C0 2 - 0.05 N^ 

60° SPHERE -CONE 

U* = 10.0 km/s 

°° 3 

= 0.01 kg/m 

T* = 200 K 


T* = 200 K 
w 



BLSSW (REF. 30) 

BL (REF. 30) 


l 


strip integral method which accounts for the upstream influence 
of the sonic corner. The present formulation of the VSL equa- 
tions cannot account for this influence. The surface pressure 
distributions corresponding to the two methods of analysis are 
shown in Figure 6. The influence of the sonic corner is sig- 
nificant only in the region 1 . 6 £ r < 2. For r<_ 1.6, the maxi- 
mum difference between the two methods of analysis is less 
than 4%. 

Convective heat-transfer rate distributions corresponding to 
the integral BL and BLSSW solutions and the present VSL solu- 
tion are presented in Figure 7. The present solution corre- 
sponds to assuming instantaneous transition from laminar to 
turbulent flow with either Equation (23) or Equation (25) used 
for the BL edge definition and Equation (22) for the normal 
intermittency factor. The convective heat-transfer rate cor- 
relation formula used in the BL and -BLSSW solutions includes a 
transition correction. The BLSSW and VSL solutions show the 
expected influence of vorticity interaction, but the VSL solu- 
tion predicts heating rates that are approximately 30 % lower 
than that given by the BLSSW solution. These differences are 
approximately the same as those obtained in the region of strong 
vorticity interaction for the previously discussed perfect gas 
solutions. For the present problem, these differences are not 
excessive and are generally to be expected when comparing dif- 
ferent numerical techniques using different turbulence models. 
For the laminar flow region, the different solution procedures 
are in excellent agreement. 

VSL solution results are presented in Figures 8-11 for a 40° 
half-angle hyperboloid for freestream conditions corresponding 
to a typical Jovian entry trajectory point. Wall shear-stress 
and convective heating-rate distributions are compared with 
the solution obtained using the BL analysis of Bartlett and 


27 



NONRADIATING EQUILIBRIUM GAS (0.74 H 0 - 0.26 H ) 

2 e 

U* = 39.09 km/s 
P* = 3.B9 xio" 4 kg/m 3 



Figure 8. 'Wall shear stress distribution for a 40° hyperboloid. 



NONRADIATING EQUILIBRIUM GAS (0.74 H 0 - 0.26 H ) 

2 e 

U* = 39.09 km/s T* = 4564 K 

°° -4,3 w 

p* = 3.89x10 kg/m r* = 0.229 m 

00 n r 

T* = 145 K N =7.4xlO P 

°o Re 

200 


CONVECTIVE 
HEATING ioo 

-qJ^TE, 2 

q c, w, MW/m 

0 0.3 1.0 1.3 2.0 2.5 3.0 

s 

Figure 9. Heating-rate distribution for a 40° hyperboloid. 



ro 

^0 





to 

H 


Figure 11. Temperature profiles at s = 2 for a 40° hyperboloid. 
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Kendall. The inyiscid constant entropy edge conditions for 
this BL solution were determined using the TFD procedure of 
Reference 24. 

Wall shear-stress and convective heating-rate distributions 
are shown in Figures 8 and 9, respectively. With the exception 
of the region near the point of instantaneous transition from 
laminar to turbulent flow, the wall shear-stress distributions 
for both laminar and turbulent flow show differences of less 
than ±10%. Laminar convective heating-rate distributions are 
essentially identical and are shown as a single curve, and for 
turbulent flow, the two methods of analysis differ by ±5% 
except in the immediate vicinity of the transition point. 

VSL velocity profile data are compared with the TFD inviscid 

24 ^1 

solution and the turbulent BL solution in Figure 10. In 

the outer predominantly inviscid portion of the shock layer, 
the VSL solution is in good agreement with the inviscid solu- 
tion. The velocity profile data in the inner viscous region 
shows satisfactory agreement near the surface but large dif- 
ferences are obtained near the outer edge of the BL. The tem- 
perature profiles (see Figure 11) are in ' satisfactory agreement 
with maximum differences of less than ±10% in the viscous inner 
region and less than ±3% in the outer inviscid region. 

The VSL solution was first obtained using the transport property 

32 

data for hydrogen and helium presented by Lick and Emmons and 
Esch, et al.^ Using these transport properties, large differ- 
ences between the VSL solution and that of Reference 31 were 
noted, and the transport property data were changed to conform 
with those used in the BL solution. The heating-rate distribu- 
tions corresponding to the transport property data of Refer- 
ences 16 and 32 are presented in Figure 12 for both radiating 
and nonradiating laminar and turbulent flows. The convective 
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Figure 12. Heating-rate distributions for a 40° hyperboloid. 



heating rates obtained with these transport properties are 
approximately 25$ to 40$ greater than those presented in 
Figure 9 for the nonradiating case. These two solutions for 
identical freestream conditions and body shape effectively 
demonstrate the difficulties to be expected when attempting 
to make comparisons of different analyses for equilibrium gas 
mixtures. For this problem, radiation heat transfer from the 
shock layer results in a reduction of temperature gradients 
in the viscous inner region with a corresponding reduction in 
convective heating rates at the surface for both laminar and 
turbulent flows. The radiation heat transfer to the surface 
is not significantly influenced by turbulence for this case. 
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3. DEVELOPMENT OF AN ADVANCED TURBULENCE MODEL 

Results of the preceding section show that, using the mixing- 
length method, satisfactory solutions for viscous shock 
layers can be obtained provided the necessary boundary- 
layer parameters — boundary-layer thickness, displacement 
thickness, and edge conditions — are accurately defined. These 
boundary-layer parameters are needed for the mixing-length 
model, and are otherwise not essential to the solution of the 
VSL equations. As shown in Subsection 2.2, for nonradiating 
perfect gases and mixtures of perfect gases in chemical equil- 
ibrium, these parameters can be satisfactorily defined by assum- 
ing that the edge of the predominantly viscous layer is located 
at the point where the local total enthalpy is 99*5$ of the 
freestream value. However, since energy is transferred out of 
the shock layer by radiation, this definition is unsuitable 
for radiating gas flows. The composite dif fusion-conduction- 
dissipation index defined by Equation (25) provides a more 
suitable method for defining boundary-layer thickness, although 
the definition is adversely influenced by the nonvanishing nor- 
mal gradient of the tangential velocity component in the nose 
region of blunt bodies. Hence, a turbulence modeling technique 
which is independent of the usual boundary-layer parameters is 
needed for VSL applications. 

Recent progress in development of phenomenological turbulence- 
model equations for boundary layer applications indicate that 
these techniques can be appropriately modified for VSL appli- 
cations. The most extensively developed and tested turbulence 

OO Oil 

model of this class is that proposed by Saffman 5 and 
Wilcox. 35»36 e q Ua tions have been formulated for compres- 

sible flows, ^“36 an q have been shown to be quite accurate 

in predicting f reestream-turbulence-induced stagnation heating 
37 

augmentation. To achieve a solution, the model requires no 
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advance knowledge of a given turbulent flow aside from boundary 
and initial conditions. This fact, combined with the model's 
independence of boundary- layer parameters, illustrates the model's 
suitability for VSL applications. 

In this study, preliminary steps have been taken to incorporate 

O f. 

the model equations postulated by Wilcox and Traci J into the VSL 
code described in Section 2. The model equations are stated in 
Subsection 3.1- Then, in Subsection 3*2, boundary conditions 
suitable for surfaces with mass injection are developed. Model- 
predicted flow structure near the edge of the predominantly vis- 
cous region is then analyzed in Subsection 3*3- Finally, Sub- 
section 3*4 proposes an algorithm for reducing computing times 
attending solution of the model equations. 


3.1 THE MODEL EQUATIONS 

The mean equations of motion. Equations (1) through (6), remain 
unaltered. The eddy viscosity, y T , is defined in terms of the 
turbulent mixing energy, e, and the turbulent dissipation rate, 
w, as follows: 

y T = p 2 e/w (30) 

In the Wilcox-Traci model, the mixing energy and dissipation rate 
satisfy the following nonlinear diffusion equations. 


u 3e . 3e 
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[y(l + a*e+)ff] ♦ } (31) 


p(l^!f 2+ 0 - {-PS - [ e + 2a(|A) 2 _ 


w > w‘ 


+ a 




j cos e 


1+nK r+n cos 0 
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( 33 ) 


The quantities S and & are given by the following: 

i = pe^/w , s = /sTTsTT 

ij Ji 


where S. . is the mean rate of strain tensor. There are six 

closure coefficients appearing in Equations (31) and (32), 

viz, a,a*', B , 3* ,c , and a*. The following values have been 

established by arguments based on general properties of 

36 

turbulent flows . 


M 20 p 100 
» ■ I ? 

a = ~[l - — exp (-e+/2)] 

a * = ]^[ 1 - II exp (-2e + ) ] 


(34) 


3.2 SURFACE BOUNDARY CONDITIONS 

34 

At solid boundaries with no mass addition, Saffman and Wilcox 
have shown that in addition to the no-slip velocity and surface- 
temperature (or heat transfer) boundary conditions, the follow- 
ing boundary conditions must be imposed on e and w: 

p 2 u 2 

e = 0 , w = — 4s R (k + ) at n = 0 (35) 

M w °° 


where S^ is a universal function of the roughness height of 
the surface, k, and 



w 


u T k/y 


w 


The constant a* is 3/10, the high-Reynolds-number (e + >>l) 
value of a*. For boundary layers with surface mass injection 
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(blowing), the introduction of an additional velocity scale 
(v = normal flow velocity at the surface) suggests that some 
modification of Equation (35) may be required for flows with 

O O 

blowing. Andersen, et al , provide further evidence that 
the dissipation-rate boundary condition must be revised when 
blowing is present by showing, from correlation of their 
experimental data, that the law of the wall assumes the 
following modified form: 

u 1 u r y 

= ^T - lo S ~ ~ + C (36) 

u T k]_ v 

The modified "constants" k^ and C are related to k^, C, and 
the blowing velocity as follows: 

- VCl + 7.7 v+) 

C = C-50vS 7 2 ( v* ) 2 

Since Saffman and Wilcox’s computations show that C is strongly 
affected by the value of S D , Andersen's data indicate that a 

n 

modification to the model equation boundary condition is 
required for computing blown boundary layers. Since mass 
injection is of key importance in VSL applications, we thus 
consider effects of blowing on model-predicted flow structure. 



As in the case of rough surfaces, we can establish suitable 
boundary conditions for flows with mass injection by analyzing 
the viscous sublayer of an incompressible plane-wall boundary 
layer. Denoting distance normal to the surface by y, the 
blown sublayer equations are 


v(l+e + ) 


du 

dy 
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+ V u 
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(38) 


v w 


de 

dy 




+ 



l+a*e + 



(39) 
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where v is kinematic viscosity. 


(HO) 


To establish boundary conditions for this fifth-order system, we 
assume the effluent gas is free of turbulent fluctuations so that 
e vanishes at the surface. The horizontal velocity, u, vanishes by 
virtue of the no-slip condition, and in analogy to Equations (35) 
we write 

pu x . 

u = e = 0 , w = ^ S (v+) at u y/v = 0. (41) 

oo 


where Sg is a universal function of v+. Examination of wall-layer 
structure (i.e., the asymptotic solution for u T y/v- > - 00 ) predicted 
by Equations (38-40) yields two more boundary conditions valid as 
u y/v -> °°. As can be easily shown, the turbulent diffusion term in 

T 

the mixing-energy equation is negligible for u^y/v-* 00 , and there 
follows 


e 





+ 


w 



w 


p du 
a* dy 

oo u 


as 


u T y/v ->• 00 


(42) 


In addition to establishing boundary conditions, further examina- 
tion of the wall-layer solution demonstrates direct correspondence 
between model-predicted and measured effects of blowing on a tur- 
bulent boundary layer. Specifically, in the limit of weak blowing 
(v+ « 1), expanding in powers of v* shows that the velocity obeys 
a modified law of the wall similar to Equation (36). The effec- 
tive Karman constant is predicted to be 



kj/d + 


5 v + ) 
w' 


(43) 


where 


4k C-l 

2 = -^~ + 



log 


V 
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(44) 
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The functional dependence in Equation (43) is similar to that 
quoted by Andersen [Equation ( 37 ) U - Table 1 shows the variation 
of E with u y/v for k-^ = 0.4l and C= 5-5. As shown, the predicted 
value of E is reasonably close to Andersen’s value of 7-1, 
particularly at the larger values of u^y/v. 


Table 1. Variation of E with u T y/v. 


u x y/v 
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6.23 
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6.65 


Using an implicit, time-marching numerical method, sublayer calcu- 
lations were performed for several blowing rates ranging from 
v+ = 0 to v+ = 0-7393- The value of S B was varied for each blowing 
rate until close agreement with the velocity profile data of 

o Q 

Andersen, et al, was obtained (Figure 13). As expected from the 
wall-layer analysis, slopes of the various profiles (i.e., k-^ -1 ) 
are accurately predicted. Note that this means that, similar to 
the case of surface roughness, the dissipation-rate boundary con- 
dition primarily determines the variation of C with v + . Figure 14 

w 

presents a correlation of S R with v + ; an accurate analytical 

.D W 

representation of the correlation is 


S 


B 


6 

v + (l+v + ) 
w w 


> 0 
w 


(45) 


Although effects of suction (v + <0) have not been considered 

W 39 

here, computations performed by Wilcox imply that w is 

unaffected by suction. Therefore, Equation (45) should only 

be used for v + > 0 . 

w 
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Figure 13- Comparison of computed and measured sublayer velocity profiles for 
boundary layers with surface mass injection. 
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3.3 TURBULENT-NONTURBULENT INTERFACE STRUCTURE 

We now consider model-predicted flow structure near the outer 
edge of the predominantly viscous region. In solving the 
model equations, numerical difficulties have been encountered 
near the outer edge both in conventional boundary- layer com- 
putations and in preliminary VSL applications. The difficulties 
stem from existence of a sharp turbulent-nonturbulent interface 

at the boundary-layer edge. The sharpness, first analyzed by 
33 

Saffman, manifests itself as discontinuities in 9u*/9y and 
9w/9y. Accurate numerical solutions require accounting for 
these discontinuities; therefore, we review interface structure 
and propose a tentative remedy in this section. 

3.3-1 Interface Structure for a Conventional Boundary Layer 

As can be easily verified for the special case e = w= 0 in the 
freestream, production and dissipation terms [i.e., the first 
terms in brackets on the right-hand sides of Equations (31) 
and ( 32 )J and molecular diffusion are negligible near the 
interface; approaching the interface from below, the asymptotic 
behavior of e, w, and u* is 

e ^ (6-y) 2 j 

w ^ (S-y) l as y 6 ~ (46) 

u* ^ (<5-y) 1 

where 6 is boundary layer thickness. Note that the behavior 

predicted in Equation (46) is consistent with experimental 

40 

observations. Specifically, Nee and Kovasznay have observed 
that the mean velocity varies linearly with distance from the 
interface, i.e., u* ^ (S-y). Also, the length scale has been 
observed by many experimenters to approach a nonzero value 
near the edge of a turbulent flow which implies ^ pe. 2 /w^ con- 
stant as y -*■ 6 . 
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Of course such discontinuities cannot exist in a truly viscous 
flow. As shown by Saffman' 5 ' 5 and by Wilcox and Traci, there is a 
thin region centered at the boundary layer edge in which viscous 
effects become important. This region is known as the viscous 
superlayer; its role is to remove the discontinuities in 3u*/9y 
and 9w/3y. Superlayer thickness is proportional to v*/v* where 
v* is entrainment velocity. Although v* is generally very small 
for a conventional boundary layer, resolving the superlayer is 
impractical in a numerical computation. Consequently, solutions 
to the equations appear to be dispontinuities even in a viscous 
computation. 

In prior boundary layer computations, we have used the known alge- 
braic behavior of e, w, and u* to locate the boundary-layer edge 
and to obviate numerical difficulties attending the solution 
"discontinuities"; the procedure is as follows. First, the equa- 
tions of motion are integrated to a point thought to be beyond 
y= 6. The solution is then examined to determine the value of y 
at which any one (or more) of the three quantities e, w, or u* is 
less than its boundary-layer-edge value. The interface discon- 
tinuities always cause such undershooting, and the dissipation 
rate is usually the first quantity to undershoot its edge value. 

If no such undershoot occurs, the equations are integrated to a 
larger value of y. Ultimately, the boundary-layer edge is located 
by linear extrapolation, from points below the undershoot loca- 
tion, of the computed values of either u* or w. For classical 
boundary layers, this procedure obviates numerical difficulties 
attending insufficient resolution of the superlayer; however, the 
procedure’s applicability for YSL's is uncertain. In fact, a key 
advantage of using the model equations, i.e., independence of 
boundary layer thicknesses, is lost if this procedure Is used. In 
the following sections we analyze interface structure for stagna- 
tion point flow and postulate an alternative procedure. 
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3-3-2 Interface Structure for Stagnation Point Plow 


To determine similarities between interface structure for VSL's 

and conventional boundary layers , we now analyze flow near a 

stagnation point (Figure 15). This flow is relevant because a 

VSL begins at a stagnation point and approaches a conventional 

boundary layer farther downstream. For perfect-gas flow near 

the stagnation point of a bluff body, the equations of motion 

37 

admit a similarity solution; the similarity variables are 


5 = 



q 


V 


T* 1 


y 


(47) 


The quantity c* is the inviscid flow velocity gradient which 
depends on the shape of the body and the Mach number. Provided 
the shock standoff distance is large compared to the thickness 
of the viscous region, the mean flow outside the viscous region 
(y>6) satisfies the usual inviscid solution near a stagnation 
point given by 


u* = c*x , v* = -(j+l)c*y 


(48) 


where j=0 for two-dimensional flow and j =1 for axisymmetric flow. 
The various dependent variables are scaled as follows: 


u* = 

v*c* 5 U(n) 

II 

* 

> 

v*c* v(n) 

T* = 

T «»T(n ) 

e = 

v*c* E(n) 

w = 

p * c * w ( n ) 

itions 

of motion are 


(49) 


I 
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Mass Conservation 


g + (J+Do - o 


(50) 


Horizontal Momentum Conservation 


V— = 1 - U 2 + — [~(l+e + ) 

dn dn L dn J 


(51) 


Energy Conservation 


dT 

dn 


_d_ 

dn 


( 


N Pr N Pr,T 


dT 

dn 


Turbulent Mixing Energy 

V.f| = {2a*/2j+T | U | - g*w} E + d_[(l+a«e + ) f|] 


(52) 


(53) 


Turbulent Dissipation Rate 

= {2a/2j4T 1 U | - BW- 2aw(|^) 2 } W 2 

where L and e + are defined as 



[< 1+ « + >ff] <54> 


L = eVw ) 
e + = E/W ) 


(55) 


Equations (50)-(55) constitute a 9th order set of coupled 
ordinary differential equations. For a perfectly smooth wall 
with no mass injection, the boundary conditions are 

U(0) = E ( 0 ) = V(0) = 0 \ 

T( 0 ) = T*/T* ( (56) 

w(n) ■+ n" 2 as n o ) 
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while (denoting Re^ = p*U*D*/jj* where D* is body diameter) 


u(n) 

» l , T(n) + i , 

1 


E(n) 

- Re D (e„/nJ)(uycD) 

> as n -* 00 

(57) 

w(n) 

- ' 5 ? <»='•»> 1 

CO 

1 



Unlike the conventional boundary layer, the production rates in 
Equations (53) and (54), i.e., 2ot*/2j +1 | U | and 2a/2j +1 ' | U | , are 
nonvanishing above the viscous region. As a consequence, there is 
apparently no sharpness at the turbulent-nonturbulent interface. 
For example. Figures 16-18 show computed e, £, and e + profiles 
(i.e., the curves labeled "without intermittency " ) for a case in 
which both e and w are very small outside the viscous region 
(y//v*/c*‘ ^ 6 ); computations have been performed with the same 
implicit, time-marching method used in the sublayer computations 
of Subsection 3-2. All three quantities smoothly approach their 
freestream values. Thus, the numerical difficulties present in 
our early VSL computations are not associated with flow in the 
immediate vicinity of the stagnation point. Rather, these diffi- 
culties must be caused by appearance of the sharp interface as the 
YSL approaches conventional boundary-layer structure. 

3 - 3-3 Tentative Modification to the Model Equations 

Having identified the source of difficulty as the sharp interface, 
we now propose a tentative remedy. Inspection of Equations ( 31 ) 
and ( 32 ) shows that, unlike the basic conservation equations, we 
can never obtain solutions for which e and w simultaneously 
approach constant, nonvanishing values. To understand this, note 
that vanishing of the source term in either equation precludes 
vanishing of the other equation’s source term. As shown by Wilcox 
and Traci,' 3 , the source terms play a secondary role right at the 
interface. Nevertheless, the stagnation-point analysis above 
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Figure 16 . Computed variation of turbulent mixing energy near a 
stagnation point. 

J=r 

VO 











shows that alterations to the source terras give rise to signifi- 
cant modification to interface' structure. 


As a numerical experiment and a potential remedy for the numeri- 
cal difficulties encountered in YSL computations with the model 
equations, we propose forcing the source terms in Equations (31) 
and (32) to vanish outside the predominantly viscous region. To 
accomplish this end, we in effect introduce intermittency factors, 

I and I , which multiply the source terms, i.e. 

0 w 


[a*pS - £*wje -* I e £a#pS - g*wje 

Japs - [e + japs - [. + 2p(H) 2 ]>'j 


(58) 

(59) 


A VSL computation was performed with 


I 


e 


I = t/t 
w w 


( 60 ) 


where r is Reynolds shear stress. No numerical difficulties 
were encountered, thus verifying that the nonvanishing of the 
source terms above the viscous layer is the primary source of 
difficulty . 

In order to achieve a more universal (tensor invariant) formula- 
tion, an alternate definition for I and I must be proposed. 

One possibility follows from noting that 8£/9n vanishes above 

both conventional boundary layers and above the 

viscous region in VSL’s. Thus, we propose the following 


I 


e 


I = 
w 



(61) 


To assess the effect of the intermittency modifications on 
predicted stagnation-point flow properties, the computations 
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discussed in Subsection 3.2.2 were repeated using the intermit- 
tency defined in Equation (6l) with r = 50 • Figures 16-18 show 
that the solution is strongly affected by int ermittency throughout 
the flow, including deep in the viscous region. Because the 
eddy viscosity is so small, the mean velocity profile is 
unaffected. However, the intermittency modification conceiv- 
ably has a nontrivial effect on velocity profiles in regions 
where the flow is fully turbulent. Additionally, the model’s 
ability to predict transition may be affected by damping the 
source terms. Hence, the intermittency modification must be 
viewed as tentative until the impact of the modification upon 
these and other important aspects of model-predicted flow 
properties have been quantified. 

3.4 COMPUTING-TIME REDUCTION 

Even with the intermittency modification, computing times with 
the VSL program are excessive when the model equations are 
used. The maximum streamwise step size at which converged 
solutions can be obtained is less than a tenth of that which 
can be used in mixing-length computations. In order to reduce 
computing time, a special integration algorithm has been 
devised. In this section, we discuss the reason for the 
reduced step size and present the integration algorithm. 

Step size must be reduced with the model equations because the 
source terms are very large and negative close to a solid bound- 
ary. Consequently, in a boundary-layer marching procedure, there 
is a tendency to overestimate the decay of numerical errors from 
one station to the next (Figure 19). Unless very small steps are 
taken as we march from x = x-^ to x = X £ » the magnitude of a numeri- 
cal error, fie, can grow very rapidly and ultimately can destroy 
the solution. By taking smaller steps, we more nearly approach 
the rapid exponential decay implied by the large and negative 
source terms, thus preserving solution accuracy. 
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Above the viscous sublayer the source terms are generally much 
smaller. Hence, it may be necessary to take small steps in 
the streamwise direction only in the sublayer while much 
larger steps may be permissible elsewhere. To test this pos- 
sibility, two incompressible flat-plate boundary layer compu- 
tations were made with a second-order-accurate, implicit 

boundary-layer code. In the computations, we solved the con- 

t h 

servation equations only at every n streamwise station. At 
each n th station, the conservation equations were solved all 
the way across the boundary layer. Below a specified distance 
normal to the surface, y-? , the model equations were solved 

at every streamwise station. In addition, the model equations 

were solved all the way across the boundary layer at each 

th . , . 

n station. 

The first round of computations was for a fully turbulent 

boundary layer; the computation used 100 mesh points between 

y = 0 and y = 6 while the value of n was 8. Figure 20 shows the 

ratio of computing times, t/t , for three values of j : t is 

o iricLx. 

the computing time required when n = 8 and t Q is the computing 

time required for a normal boundary layer computation, i.e., 

n= 1 and j = 100. As shown, there is a 17$ reduction in 

computing time when the model equations are solved up to y = 6, 

i.e., for j = 100. However, when j Is 30, computing time 
max max 

is nearly halved. Further reduction of J max results in a 
slight increase in computing time; this increase is caused by 
an attendant increase in the number of iterations required to 
obtain a converged solution. In all computations we observed 
no deterioration in solution accuracy. 

The second round of computations was for a laminar/transitional 
boundary layer; starting from the plate leading edge, the equa- 
tions of motion were integrated up to the model-predicted 
transition point at a plate-length Reynolds number. Re , of 

X 
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2.1*10 6 . Using j =100, we varied n from 1 to 25. As shown 

max 

in Figure 20, t/t drops by 30$ as n increases to 10; a slower 
decrease in t/t follows as n increases further with a large n 
asymptote of t/t Q =0.60. Varying n causes no variation in 
transition Reynolds number. 

This new integration algorithm, in effect, permits an order of 
magnitude increase in streamwise step size. The reduction in 
computing time is .somewhat disappointing however. A signifi- 
cant increase in the number of iterations required to achieve 
converged solutions prevents computing time from decreasing 

proportionally with j and/or n. Nevertheless, the compu- 

max 

tations verify the proposed cause of step size reduction, viz, 
the large and negative source terms in the model equations, 
particularly near a solid boundary. Having attained such 
understanding, we are now in a position to develop even more 
efficient integration algorithms. 
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4. CONCLUDING REMARKS 


As shown in Section 2, numerical solutions for the viscous-shock- 
layer equations compare favorably with inviscid flowfield solu- 
tions in the inviscid outer portion of the shock layer. Favorable 
comparisons with boundary-layer solutions are obtained for the 
predominantly viscous inner portion of the shock layer both near 
the nose and far downstream where vorticity interaction is negli- 
gible. In flow regions*where vorticity interaction is strong, the 
viscous-shock-layer solutions appear to be more satisfactory. 

First-order boundary-layer theory with corrections for streamline 
swallowing appear to overpredict boundary-layer and displacement 
thicknesses. For turbulent flow, this results in an apparent 
overprediction of skin friction and heat transfer if vorticity 
interaction is strong. This behavior is the result of the implied 
nonvanishing normal gradients at the edge of the boundary layer 
when longitudinal entropy gradients are considered. 

Experimental data are needed for hypersonic turbulent flows over 
blunt bodies to assess the accuracy of the present numerical solu- 
tion procedure and to determine the validity of mixing-length 
turbulence models in the presence of strong vorticity interaction. 

Turbulence models that are independent of conventional boundary- 
layer thickness parameters are needed for both the present method 
of analysis and for numerical solution of the Navier-Stokes equa- 
tions. The analyses of Section 3 demonstrate the potential of two- 
equation turbulence models in the YSL context. While further 
development of the model and methods for efficiently integrating 
the model equations will be needed to achieve a practical engineer- 
ing design tool, the results are very encouraging. 
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